//
// Created by tview on 20-12-22.
//

public struct complex//定义复数
{
    public double real;
    public double imag;

    // 写个函数用于显示
    public void show()
    {
        if (Math.Abs(real) < 0.0001) real = 0;
        if (Math.Abs(imag) < 0.0001) imag = 0;
        if (imag > 0) Debug.Log(string.Format("{0} +{1}i", real, imag));
        else if (imag == 0) Debug.Log(string.Format("{0}", real));
        else Debug.Log(string.Format("{0} {1}i", real, imag));
    }
}

public static complex[] calDFT(double[] f)  //(信号，信号长度)
{
    int N = f.Length;
    complex[] F = new complex[N];
    //离散不傅立叶变换，其实达到最大取值点，就完全满足京都了。时间复杂度是ｏ(ｎ)。快速傅立叶变换未0(n*log(n))
    for (int n = 0; n < N; n++)
    {
        // 声明
        F[n].real = 0;
        F[n].imag = 0;
        for (int t = 0; t < N; t++)
        {
            // 计算 exp(-i * 2PI * n / N * t)
            complex temp;
            // 欧拉公式 exp(ix) = cos(x) + isin(x)
            temp.real = Math.Cos(-2 * Math.PI / N * t * n) * f[t];
            temp.imag = Math.Sin(-2 * Math.PI / N * t * n) * f[t];

            F[n].real = F[n].real + temp.real;
            F[n].imag = F[n].imag + temp.imag;
        }
    }
    return F;
}

public static complex[] calIDFT(complex[] F)  //(频谱)
{
    int N = F.Length;
    complex[] f = new complex[N];
    for (int t = 0; t < N; t++)
    {
        // 声明
        f[t].real = 0;
        f[t].imag = 0;
        for (int n = 0; n < N; n++)
        {
            // 计算 exp(i * 2PI * n / N * t)
            complex temp;
            // 欧拉公式 exp(ix) = cos(x) + isin(x)
            double real = Math.Cos(2 * Math.PI * n / N * t);
            double imag = Math.Sin(2 * Math.PI * n / N * t);
            // 复数乘法
            temp.real = F[n].real * real - F[n].imag * imag;
            temp.imag = F[n].real * imag + F[n].imag * real;

            f[t].real = f[t].real + temp.real;
            f[t].imag = f[t].imag + temp.imag;
        }
        f[t].real = f[t].real / N;
        f[t].imag = f[t].imag / N;
    }
    return f;
}


